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Abstract 

We present a comparative study of several dynamical systems of increasing complex- 
ity, namely, the logistic map with additive noise, one, two and many globally-coupled 
standard maps, and the Hamiltonian Mean Field model (i.e., the classical inertial 
infinitely-ranged ferromagnetically coupled XY spin model). We emphasize the ap- 
pearance, in all of these systems, of metastable states and their ultimate crossover 
to the equilibrium state. We comment on the underlying mechanisms responsible 
for these phenomena (weak chaos) and compare common characteristics. We point 
out that this ubiquitous behavior appears to be associated to the features of the 
nonextensive generalization of the Boltzmann-Gibbs statistical mechanics. 
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1 Introduction 



The statistical-mechanical comprehension of the occurrence and eventual dis- 
appearance of out-of-equilibrium quasistationary states (QSS), referred here 
also as metastable states, is an important undertaking that nowadays attracts 
the attention of researchers in various fields, ranging from weakly chaotic non- 
linear dynamics to slow glassy dynamics in condensed matter. Under specific 
classes of initial conditions, the dynamics in these systems displays two-step 
relaxation processes, entering first one or more periods of little variation before 
a crossover to the final equilibrium state. The QSS are conjectured to emerge 
as the result of long-ranged correlations or interactions amongst the elements 
of the system or long-ranged time correlations in the orbits of iterated maps. 
An interesting statistical-mechanical model that displays QSS is the Hamil- 
tonian Mean Field (HMF) model [1], that in its simplest version consists of 
a set of N inertial XY classical spins or rotors all of which interact equally 
with each other. But also, interestingly, the most relevant features of QSS 
dynamics are displayed by considerably simpler nonlinear dynamical maps, as 
it is the case of conservative coupled maps [2] and also by single, dissipative, 
one-dimensional maps, like the prototypical logistic map close to the edge of 
chaos and in the presence of external noise [3]. 

In this paper we present a comparative study of the aforementioned models, 
starting our discussion with low-dimensional nonlinear maps and ending it 
with the HMF model. We emphasize the analogies found between their dynam- 
ical properties and indicate the mechanisms that give raise to the QSS as well 
as those responsible for the crossover to the equilibrium state. In all of these 
cases the quantities that describe dynamical evolution seem to be associated 
to the predictions of the nonextensive generalization [4] of the Boltzmann- 
Gibbs (BG) statistical mechanics. Also, in all of these examples there is clear 
evidence of the non-commutability of the infinite time limit with other impor- 
tant parameter limits, such as the thermodynamic limit (N — > oo) or the limit 
of vanishing noise amplitude. Moreover, there is a common indication that, 
in the characterization of the phase space available for dynamical evolution, 
there appears evidence of fractalization and of ergodicity failure. 

The paper is organized as follows. In the next section we describe how the 
dynamics at the edge of chaos of the logistic map leads, via the addition of 
external noise, to the crossover from QSSs to a strongly chaotic state that 
is the analog of the BG equilibrium state [3]. We present here new numerical 
evidence by means of an ensemble implementation of the dynamics in this sys- 
tem. In section 3 we note that the emergence of QSSs and eventual departure 
from them appear associated to the complex structure of phase space accessi- 
bility in symplectic maps governed by the Kolmogorov-Arnold-Moser (KAM) 
theorem [2]. We also show new results that indicate that the QSS, present at 
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the lowest possible dimensions, persists as hundreds of symplectic maps are 
(globally) coupled. Section 4 is devoted to the QSS of the HMF model. Par- 
ticular attention is paid to the time evolution of the anomalous probability 
distribution functions (PDFs) for the velocities, and also to basic questions 
referring to the measurability of temperature and to the zeroth principle of 
thermodynamics as it relates to the QSS [5]. Finally, in section 5 we draw our 
conclusions, we make reference to the existing analogies amongst the different 
models described and discuss the underlying mechanism that gives raise to 
the QSS. 



2 QSS in the logistic map with external noise 

Let us start our analysis by recalling [3] that phenomena like two-step relax- 
ation and aging can be displayed at the lowest possible dimension, namely, by 
means of a one-dimensional map. For this purpose we consider the paradig- 
matic logistic map with additive white noise, 

x T+1 = 1 - fix 2 T + £ T a, r = 0,l,..., xe[-l,l], A*e[0,2], (1) 

where £ r is Gaussian-distributed and delta-correlated, (£t)c = 0> (£t£t')£ = ^tt', 
and a measures the noise intensity. As it is well known, for a = the Feigen- 
baum attractor at the edge of chaos /i — /i c = 1.401155189... is the accu- 
mulation point of both the period-doubling and band-splitting cascades. The 
Lyapunov coefficient is zero [6] and the attractor is a cantor set of fractal 
dimension df = 0.5338.... It has been analytically proved by means of the 
conventional renormalization-group approach for this type of map that the 
dynamics at this critical point is unmistakably associated to the nonexten- 
sive statistical mechanics [7]. Specifically, the sensitivity to initial conditions 
is actually given in a closed form by a g-exponential function and it is related 
to the rate of entropy production via a g-generalized Pesin identity, linking 
the sensitivity to initial conditions to the nonextensive expression for the en- 
tropy, S q = ln q W (where W is the number of equally-probable positions of 
an ensemble of iterates at a given time). For this stationary state the in- 
dex q takes the specific value q = 0.2445.... Let us remind the reader that 
the g-exponential and its inverse the g-logarithm are defined, respectively, as 
expire) = [1 + (1 — q)x] 1 / 1 ~ q and \n q (x) = (rr 1 " 9 — 1)/(1 — q), and that these 
expressions reduce to the usual exponential and logarithmic functions in the 
limit q — > 1. 

More precisely, the Feigenbaum attractor can be described by means of a set 
of position subsequences generated when the iterative map is given the initial 
position xo = [7]. See the empty circles in Fig. 1(a) where we have plotted 
the absolute value of these positions \x T \ and time r at which they are reached 
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Fig. 1. (a) Phase space structure of the logistic map of the logistic map with 
and without additive noise. Empty circles correspond to the attractor at n = fi c 
and a = 0, obtained by iterating xo = (the numbers label time r = 1, 16). 
Small dots correspond to xo = 0.56023, close to a repeller, the unstable solution 
of x = 1 — /i c x 2 . Full (and dashed) lines represent two trajectories for fi = fj, c , 
a = 10~ 3 (and a = 10 -6 ). (b) Plateaux structure. Time evolution (using the shifted 
time t = 2 n — 1) of the variance of an ensemble of M = 10 4 trajectories all starting 
at xo = 0, for different orders of magnitude of the noise strength a. Notice that 
the crossover times coincide with o r ~ l (with r ~ 0.633) indicated with arrows (see 
text). 



in logarithmic scales. Each of the subsequences is labelled by an integer k = 
0, 1, ... and each is generated by the time subsequence = (2k + l)2 n ~ k , with 
n — 0, 1, ... and n> k. This can be appreciated in Fig. 1(a). In the following we 
will focus on the principal subsequence tq = 2 n . Another important feature of 
Fig. 1(a) is that the logarithm of the attractor positions appear clearly divided 
into an infinite number of equally-spaced horizontal bands. In the empty space 
between any two bands there are points called repellers, the unstable solutions 
of the iterated equation x = 1 — fi c x 2 . (See the small dots in Fig. 1(a)). At 
the chaos threshold of the logistic map with a ^ the dynamics is governed 
by two competing phenomena: On the one hand the action of the noise is to 
drive the iterates out of the bands, and on the other hand the instability in 
the neighborhood of the repellers sends the positions back into the attractor 
bands. When the noise intensity is large enough to fill a gap between bands, 
the result is that the iterates are able to jump from one band to the next. 
(See full and dashed lines in Fig. 1(a)). Since the band gaps in |x|-space are 
of variable size, increasing in size as \x\ — * 1, the motion becomes confined 
among some few bands for a certain amount of time but at longer times the 
iterates are able to spread to other more distant bands. 
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The analysis of the crossover time t c for the broadening of the attractor bands 
produced by the noise can be made quantitative. At the edge of chaos /i = fi c , 
for a 7^ 0, the iterated position of x = at shifted time t = r — 1 is given 
by (see [3] for details): 

x t = exp 2 _ (? (-A (? £)[l + &(7exp r (A r £)], (2) 

where q = 1 — In 2/ In a = 0.2445..., \ q = In a/ In 2 = 1.3236... is the general- 
ized Lyapunov coefficient, r = 1 — In 2/ In k ~ 0.633, A r = In k/ In 2 ~ 2.727 
(a = 2.50290... being one of Feigenbaum's universal constants and k ~ 6.619 
a constant associated to the noise term). Let us consider now a different tra- 
jectory x' t also starting at the origin at t = 0, x' Q = 0. Since the noise acts 
independently on each trajectory, after time t the orbit's separation is given 
by * 

x t - x\ = exp 2 _ q (-\ q t) exp r (A r ty(& - Q. (3) 

At time t = 2 n — 1 the separation between bands at which the iterate is 
located is of order exp 2 _ 9 (— X q t). Using £ t — £' t ~ 1 we obtain an estimate of 
the crossover time t c at which the noise- assisted jumps into different bands 
starts to take place. This is 

exp r (A r i c )o- ~ 1 =>• t c ~ o- r_1 . (4) 

In order to display the many-plateaux structure of time evolution in this dy- 
namical model we consider next an ensemble of M = 10 3 independent copies 
of the map (1) at the edge of chaos ji = /i c . The initial positions are xq = VM 
and we measure the variance — (x t ) 2 (which plays a role analogous to the 
temperature in the statistical-mechanical description of Hamiltonian systems) 
of the distribution of positions at time t — 2 n — 1. Fig. 1(b) shows a numerical 
corroboration of previous results [3], for different values of the noise intensity 
a. While in the first plateau the variance is zero, as the dynamics of all M 
trajectories is confined to the attractor, at the crossover time t c ~ cr r_1 the 
ensemble of orbits starts spreading and the variance increases until a second 
plateau is reached where it saturates as a consequence of the confinement in- 
duced by the unstable repeller regions. As it can be well appreciated in Fig. 
2, for long enough times it is possible to see another increment of the variance 
as the iterates are able to reach other bands of the attractor. The whole pro- 
cess ends when all bands become occupied and the variance is consequently 
of order 1. 

A two-step process of relaxation is known to be one of the main dynamical 
properties displayed by supercooled liquids close to glass formation [8]. Two 
other characteristic features of glassy dynamics are the so-called Adam-Gibbs 
formula and the scaling property known as aging. The first is a relation be- 
tween the relaxation time t c (e.g. the viscosity or the inverse of the difussivity) 
and the entropy S con f associated to the number of possible molecular minimum 
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Fig. 2. Plateaux structure for the logistic map with additive noise. Time evolution 
(using t = 2 n — 1) of the variance of an ensemble of M = 10 4 trajectories all starting 
at xo = 0, for different values the noise strength a. In order to display data with 
different orders of magnitude, the y-axis has been plotted in a logarithmic scale 
with a shift of 10 6 . 

energy configurations in the fluid [8]. The second is the loss of time translation 
invariance, named aging [9], that is due to the fact that properties of glasses 
depend on the procedure by which they are obtained. The time decrease of 
relaxation functions and correlations exhibit a scaling dependence on the ra- 
tio t/t w where t w is a waiting time. As the counterpart to the Adam-Gibbs 
formula, it has recently been shown [3] that the logistic map at fj, c for a ^ 
presents a relationship between the plateau duration t c , and the entropy S con f 
for the state that comprises the largest number of (iterate positions) bands 
allowed by noise intensity. This entropy is obtained from the probability of 
chaotic band occupancy at position x [3]. Also, the trajectories at \x c in the 
limit o — > have been shown [3] to obey a scaling property characteristic of 
aging, of the form x t+tw = g^ tw \0) exp q (—\ q t/t w ) where g( tw \0) is the t^-times 
composition of the Feigenbaum fixed-point map function g(0). 



3 QSS in coupled symplectic maps 

Symplectic maps are important tools in the study of Hamiltonian systems. If 
we consider a time-independent Hamiltonian system with n degrees of free- 
dom (2n-dimensional Gibbs T-space), a (2n — 2)-dimensional symplectic map 
is the result of taking a Poincare section over the constant energy hypersurface 
[10]. The recurrence time is discrete and the map is useful in displaying the 
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Fig. 3. QSS in the standard map (i.e., N = 1). (a) T(t) for typical values of a. 
We start with "water bag" initial conditions (M = 2500 points in < 6i < 1, 
pi = 0.5±5 10~ 4 ). Inset: Inverse crossover time t c (inflection point between the QSS 
and the BG regimes with the time axis plotted in logarithmic scale) vs. 1 / (ai — a c ) 2 ' 7 . 
(b) Time evolution of the ensemble in (a) for a\ = 1.1 (first row) and PDF of its 
angular momentum (second row), t = 0: Initial conditions; t = t\ = 500: The 
ensemble is mostly restricted by cantori; t = t2 = 10 5 : The ensemble is confined by 
KAM-tori. See [2] for further details. 

stationary, recurrent properties of the original Hamiltonian system. Interest- 
ingly enough, a 2n-dimensional symplectic map is also the result of a Poincare 
section in the phase space of a time- dependent Hamiltonian system with n 
degrees of freedom (2n + 1 dimensions), with a periodic dependence on time 
[10]. As a consequence of the symplectic structure, the Lyapunov spectra in 
the ci-dimensional phase space of the map (d being an even number) is charac- 
terized by d/2 pairs of Lyapunov coefficients, where each element of the pair 
is the negative of the other. 

A prototypical 2-dimensional symplectic map, that displays the KAM mech- 
anism for the transition from regularity to chaos, is the Chirikov- Taylor stan- 
dard map (see, e.g., [10] and references therein). Since fundamental dynamical 
processes, like for example Arnold diffusion, occur only for d > 2, it is interest- 
ing to consider more generally the case of N symplectic and globally coupled 
standard maps (d = 2JV), described, for example, by the equations: 

N 

6i{t + 1) = Pi (t + 1) + 9i(t) + £ Pm(t + 1) (mod 1), 

Pi (t + 1) = Pi (t) + £ sm[2n9 t (t)] (mod 1), 

where t = 1,2,..., % = 1,2, ...,iV, b e M is the coupling constant (in the 
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Fig. 4. QSS in two coupled standard maps (i.e., N = 2). (a) T(t) for b = 2 and 
typical values of a. We start with "water bag" initial conditions (M = 1296 points 
with < Oi, 62 < 1, and p±,P2 = 0.25 ± 5 10~ 3 ). Inset: Inverse crossover time t c vs. 
1/a 5 ' 2 . (b) Time evolution of the fractal dimension of a single initial ensemble in 
the same setup of (a). See [2] for further details. 



following we will fix b = 2), and a e K. is a parameter that we will use 
to control chaoticity; Q{ may be regarded as an angular variable and Pi as an 
angular momentum (denned on a compact set though). Notice that in order to 
describe a system with a phase space of finite size we are considering, as usual, 
the map on the torus (mod 1). This map is typically chaotic for large values of 
\a\, while for small \a\ ^0a complex phase-space structure develops of regular 
regions surrounded by a stochastic web. The stochastic web is connected for 
d > 2 (i.e., N > 1) and disconnected for the marginal case d — 2 (i.e., N = 1) 
(see, e.g., [11]). Intentionally we introduced a global coupling among the N 
standard maps: to obtain time evolution properties similar to those of the 
long-ranged interaction many-body Hamiltonian that we shall discuss below. 

Our coupled maps display a two-plateau structure with striking analogies with 
both the one-dimensional map discussed in Section 2 and with the HMF model 
that we analyze in section 4. As in the previous section, we start by considering 
an ensemble of M independent copies of map (5). In trying to deduce statistical 
properties directly from the dynamical behavior of a Hamiltonian system - 
with diagonal kinetic matrix and zero average momentum - it is common 
practice to identify the temperature with the average squared momentum per 
particle (see, e.g., [12]). Because we study here situations with nonzero "bulk" 
motion, a natural measure of temperature is given by the variance of the total 
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Fig. 5. QSS in N = 100 coupled standard maps, (a) T(t) for 6 = 2 and typical values 
of a. We start with "water bag" initial conditions (M = 500 points with < 6i < 1 
and pi = 0.25 ± 5 10~ 3 Vz = 1,2, ...TV), (b) Inverse crossover time t c vs. 1/a 7 - 5 . 

angular momentum, 



T(t) 



1 N 

- V 

N f~i 



p 2 m - (Mt)) 2 ) , 



(6) 



where () means ensemble average. As the classical microcanonical ensemble is 
based on the equal a priori probability postulate, we will call BG temperature 
the temperature associated to the uniform distribution in phase space: 



T, 



BG 



1 N 

-Y 



dpi p. 



dpi p, 



1/12 - 0.083 (ViV). (7) 



Of course, if the system is sufficiently chaotic or if the initial ensemble is suffi- 
ciently close to equilibrium (i.e., to the uniform distribution) the temperature 
(6) rapidly relaxes to Tbg as a function of the iteration time. Nevertheless, if 
the phase space presents complex structures and the initial ensemble is very 
far-from-equilibrium, partial barriers (see, e.g., [13]) may confine the ensemble 
in some limited region for quite long times, before allowing for a larger occupa- 
tion of phase space. Consequently, it appears a plateau with T = Tqss ^ Tbg- 

In close analogy with the case of section 4, for our out-of-equilibrium initial 
conditions we set all momenta Pi randomly distributed inside a properly chosen 
small interval of phase space (the so-called "water bag" initial conditions), 
while the angles 0j are distributed at random inside the whole interval [0, 1]. 
Fig. 3(a) displays the case of a single standard map (N = 1). For this marginal 
case the formation of total barriers for a < a c = 0.971635406... prevent the 
formation of a second plateau below this critical point. We also notice that 
unless a is large enough, the temperature associated to the second plateau is 




Fig. 6. Scaling of the crossover time with the number of coupled standard maps. 
Initial data are defined as in Fig. 5. (a) T(t) for a = 0.4, 6 = 2 and large values of 
N. Inset: Inverse crossover time t c vs. 1/N 1 - 2 . (b) T(t) for a = 0.4, b = 2 and small 
values of N. 



smaller than T BG . However, the 2-dimensional case is very useful for visualizing 
the mechanism for plateau formation (see Fig. 3(b)). Fig. 4(a) indicates that 
in d — 4 the topology of phase space differs from d = 2. Now total barriers for 
diffusive processes do not exist, and it is possible to have out-of-equilibrium 
initial data that leads to a dynamics with a Tqq S temperature two-plateau 
of long duration before a crossover to T BG . Interestingly enough, during the 
first plateau the geometry of the ensemble is characterized by a nontrivial 
fractal dimension (dj ~ 2.7 < 4 = d) that at equilibrium crosses over to the 
dimension of the whole phase space (see Fig. 4(b)). As shown by the new 
results presented in Fig. 5(a) for the 200-dimensional phase space obtained by 
coupling N = 100 standard maps, the familiar two-plateau structure persists 
for large values of N. Preliminary results suggest that even for large values 
of AT a first plateau develops with df/d < 1. This fact strongly reinforces the 
picture that there are nonzero-measure sets of initial conditions such that the 
system is maintained dynamically in a subset of the whole phase space, and 
only later occupies the entire allowed space, and adopts a BG-like statistics. 
This behavior is the more pronounced the closer to zero is the largest Lyapunov 
exponent. 

It is important to stress that in all these cases the crossover time to BG 
statistics diverges as a — > (or as a — > a c for d — 2), as it is shown in Fig. 
5(b) and in the inset of Fig. 4(a) and 3(a). The consequence of this is that 
if the t — > oo limit is taken after the a — > (or a — > a c for d = 2) limit, the 
temperature associated to these ensembles is forever Tqss ^ Tbg- The same 
feature is revealed if, fixing a, we analyze the scaling of the crossover time for 
a large N (see Fig. 6(a) where we have taken a = 0.4). Moreover, the exponent 
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in the scaling relation t c ~ A^ L2 is very close to what have been found for the 
QSS of the HMF model: t c ~ N [14]. It is also interesting to notice that, for 
the case we are discussing, this scaling is effective only for large value of iV 
(say, iV > 20). For small N a more intricate, richer behavior appear (see Fig. 
6(b)). In particular, for N = 5 we observe for the first time in the case of 
globally coupled symplectic maps evidence of many-plateau structures (like 
those in section 2). 



4 QSS in the Hamiltonian Mean Field model 

Along the previous lines of thought we can draw our attention now to a more 
complex dynamical scenario, that is, Hamiltonian many-body dynamics. In 
analogy with the globally-coupled symplectic maps of the previous section, 
we will here focus on long-range Hamiltonians. Long-ranged interacting sys- 
tems constitute nowadays an important subject in many areas of physics, 
e.g. astrophysics, hydrodynamics, nuclear physics, Bose-Einstein condensates, 
atomic clusters, plasma physics, and several others [15]. If interaction terms 
decay as r~ a and d is the dimension of the system, we say that a classical 
system has long-range interactions when a/d > 1. Long range interactions in- 
duce long range correlations, and, as a consequence of this, it has been found 
that these systems show, in a variety of situations, sensible deviations from 
the BG equilibrium behavior. 

In this section we will overview some known aspects and present as well new 
results for the HMF model [1]. This model has been largely considered in the 
literature, as it has been found to display a surprisingly rich behavior (see, e.g., 
[16] for a review). Here we are particularly interested in the appearance, under 
specific initial conditions (see below), of QSSs that after a certain amount of 
time cross over to the BG equilibrium. During the QSSs a variety of anomalous 
behaviors has been detected. Among them, vanishing Lyapunov spectrum, 
Levy walks and anomalous diffusion, anomalous distribution of the momenta, 
non-commutability of the infinite time limit with the thermodynamical limit 
(see, e.g., [17] and references therein), and, more recently, aging [18,19] and 
glassy behavior [20]. As some connections have already been established with 
nonextensive statistical mechanics [17,18,19], these QSSs are candidates for a 
situation that is correctly described by this formalism. 

The HMF Hamiltonian has the following form: 

N v 2 1 N 

H = K + V = ^ + —j: [1 - 008(6,-6,)}, (8) 

i — 1 i,j — 1 

where 6, G [0, 2tt) is the ith angle and pi G K. is the conjugate variable repre- 
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senting the angular momentum (with unit inertial momenta). It is an inertial 
version of the XY ferromagnetic spin model where the interaction terms cou- 
ple globally all the spins. Because of the presence of the kinetic part in the 
Hamiltonian, this model (in opposition with the usual Ising-like spin models) 
is naturally provided of a dynamics. Also note that it is common use (but not 
strictly necessary [21]) to divide the potential term by N in order to make the 
Hamiltonian (formally) extensive [22]. 

This system can be analytically solved within the BG canonical formalism, 
which predicts a second-order phase transition from a low-energy ferromag- 
netic phase with magnetization m = |m| ^ 0, where m = J2f=i m i/N (with 
nij = (cos $i, sin 0i)), to a high-energy phase with all spins uniformly dis- 
tributed on the interval [0, 2ir) (and consequently m — 0). The critical point 
is at specific energy u = u c = 0.75, and the complete equilibrium caloric curve 
T BG (u) can be exactly obtained [16]. On the other hand, it is possible to 
numerically integrate the Hamilton equations at fixed total energy, therefore 
allowing an analysis in the microcanonical setup. Results show nonequivalence 
between the two approaches for a range of energy densities below the critical 
point and for a certain class of initial conditions (which includes the so-called 
"water bag" initial conditions, characterized by a random distribution of the 
angular momenta around zero). Numerical simulations show in fact a first 
regime where the system stays in a non-equilibrium metastable state. It is 
possible to distinguish this QSS by observing the evolution of the dynamical 
temperature, here defined in the standard way T(t) = 2K(t)/N since the to- 
tal angular momentum is zero. After a short transient the system evolves into 
a first stage with a temperature Tqss < Tbg and it is only at longer times 
that it relaxes to the predicted equilibrium value. Very similarly to the case 
discussed in the previous section, the duration of the QSS diverges linearly as 
N goes to infinity [14]. 

Inside this framework, we present numerical results using a textbook con- 
struction (due to Gibbs in fact) of the canonical ensemble as a subsystem 
of the microcanonical one. It is a purely dynamical approach, in contrast to 
the analytical derivation mentioned above. We consider in fact a subset of M 
spins of the total N spins in the isolated (microcanonical) system, and study 
the temperature evolution of this subset, particularly during the metastable 
state. Temperature is defined as Tm (£) = 2KM(t)/M, where Km is the kinetic 
energy of the M spins. Different criteria may be used in order to choose the 
M subsystem. We show two examples, namely selecting the highest and the 
lowest temperature subsystem. Result are illustrated in Fig. 7 for u = 0.69 [5]: 
Two curves coming from above and below, respectively. In both cases Tm is 
seen to relax after some time to T/v, within the metastable state. Afterwards, 
both temperatures evolve sharing the same value and perform together the 
crossover to the predicted Tbg- This result is repeated for two system sizes, 
M = 100 (solid line) and M = 500 (dashed line), yielding the same behav- 
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Fig. 7. Temperature evolution of an isolated N = 5000 spin system (Eq. (8)) for 
u = 0.69 in grey line. We start with "water bag" initial conditions for the momenta 
(an "almost uniform" distribution inside the range —2 < pi < 2: See [5] for details) 
and magnetization m = 1 (0j = Vi). Hot (starting from above) and cold (starting 
from below) smaller subsystems are shown, both with M = 100 (solid thin line) and 
500 (dashed thick line). Inset: Magnification of the crossover between Tqss — 0.38 
and T BG - 0.467. 

ior. It is important to stress that in the case of the hot subsystems, during 
the relaxation to Tqss the temperature value Tbg is crossed with no signs of 
relaxation to it. Simulations are made using the 4th order symplectic Neri- 
Yoshida integrator [23] with energy conservation AU/U ~ 10 -4 , and "water 
bag" initial condition. This result gives a first step into the study of thermal 
meta-equilibrium of two systems, i.e., the zeroth principle of thermodynamics 
for the QSS and, furthermore, allows for the possibility of a generalized canon- 
ical treatment of the metastable state, which constitutes in fact an ongoing 
study [24]. 

We also present a new analysis of how the PDF of the momenta of the M- 
subset evolves with time. Results are given in Fig.8. We observe that the 
PDF for the M-subset essentially coincides with that of the total iV-system. 
After a transient of order t ~ 10 3 , the PDF stabilizes into an non-Maxwellian 
distribution that lasts up to the crossover time t ~ 10 6 . During this whole 
metastable stage the PDF is well fitted by a g-Maxwellian with a cut-off, 
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Fig. 8. Snapshots of the evolution of PDF for momenta in the ^-space, starting with 
a "hot water bag" (corresponding to dashed thin curve in Fig. 7). The staired his- 
togram shows the instantaneous PDF at time t = and t = W k (with k = 1,2, ...7), 
for M = 500 spins inside an N = 5000 spin system. In empty circles, the average 
of 10 3 realizations for N = 10 5 spins (all averages made during the QSS plateau). 
In dashed line, the g-exponential fitting curve with T = Tqss = 0.38 and q = 3.7. 
Finally, the curve in solid line shows the analytical Gaussian PDF. These last three 
curves are the same in every frame, and are plotted for reference. 

with temperature T = Tqss and q ~ 3.7. The value of q that we find is 
different from that calculated in [14]. This may be due to the fact that we are 
using slightly different (more ordered) initial conditions or to finite-size effects. 
Our present conjecture is that not only the value of q can change when the 
thermodynamic limit is adequately approached (first N — > oo, then M — > oo, 
and finally t — > oo), but also the location of the cut-off (in momenta) may 
go further away, towards infinity. After the crossover time the PDF leaves the 
QSS anomalous distribution and relaxes to the equilibrium quasi-Maxwellian 
with T = T BG (strictly Maxwellian, with no cutoff at all, only for N — > oo). 



5 Summary and discussion 

Let us summarize our results. We have reviewed and analyzed three quite 
different dynamical systems: First, the logistic map at the edge of chaos with 
additive noise, secondly a system composed of globally coupled standard maps, 
and finally the HMF model. In all three cases a common dynamical feature is 
found: The occurrence, for certain classes of initial conditions, of metastable 
states that after a certain time evolve into the stable equilibrium state. These 
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metastable or quasistationary states appear to be intimately associated to a 
partial, non-trivial, occupation of phase space, typically (multi)fractal (scale- 
free and perhaps network-like). This feature retains the system in a non- 
ergodic state (while held in reserve from entering into all regions of the phase 
space available to the system). A detailed description for this situation can be 
given for the simpler case of the logistic map with additive noise at the onset 
of chaos, i.e., the interplay between a multifractal attractor and a multifrac- 
tal repeller. This circumstance can also be recognized in the more complex 
system of two coupled standard maps, since our calculations clearly show 
the occurrence of quasistationary states and their crossover to final standard 
chaotic behavior. As we have seen too, the characteristic two-plateau feature 
is maintained for the case of hundreds of coupled standard maps, and so, we 
may propose the hypothesis that a fractal occupation holds even for a many- 
dimensional phase space. The similarities observed between globally coupled 
symplectic maps and the HMF model indicate that the same state of affairs 
possibly holds in the case of the QSS detected in the HMF model, a system 
for which a direct analysis of the T-space is much harder. 

The importance of these results is underscored by the fact that the QSSs do 
become permanent when conventional limits are taken in a specific order, e.g. 
the thermodynamic limit before the infinite time limit in the HMF model. 
There is a nonuniform convergence feature in the dynamical evolution that 
leads to an atypical stationary state. Non-commutability of limits occurs too 
in the case of the maps. For the coupled standard maps this is appreciated 
when the thermodynamic limit is replaced (for fixed a) by the limit of infinite 
number of coupled maps, or when the chaoticity parameter a tends (for fixed 
N) towards a specific value (a = for d > 2, and a c if d — 2); while in the 
case of the logistic map it is the limit of vanishing noise that plays an analo- 
gous role. Several other common features are found for the QSSs occurring in 
these different systems. For instance, in all of them the whole Lyapunov spec- 
trum vanishes, which seems to pinpoint the cause of the observed anomalies. 
Also, aging and other glassy properties are both present in the logistic map 
with additive noise and in the HMF model. These common features suggest 
a deeper, basic, connection among these systems. As we have seen, the usual 
BG statistical theory appears to be inadequate in explaining the metastable 
states, and that a generalization of the usual BG statistics is required. In this 
respect, the nonextensive theory [6] comes out as the strongest candidate to 
meet this test. 
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